Using waterNETN R pkg.

Using the waterNETN R package

Getting started

Installation

Step 1. Install R, RStudio, and RTools in Software Center

Step 2. Install devtools package in R:

install.packages('devtools')

Step 3. Install waterNETN from github

Note that whenever the waterNETN package is updated, you can rerun this code to install the latest version.

library(devtools)
install_github("KateMMiller/waterNETN")

Step 4. Load waterNETN R package

library(waterNETN)

Step 5. Import data

Note that R is not able to connect to files on Sharepoint or MS Teams (b/c Teams also stores all files on Sharepoint). That means you need to store data package files on your local machine or on a server (e.g. NETN Z drive). The default option for importing data will add the data package views (i.e., flatfiles) to an environment called VIEWS_WQ to your Environment work space (i.e. Environment tab in top right panel). If you would rather import each individual view into your R session, specify with the new_env argument (e.g., importData(new_env = F)).

Option 1. Import data via .csv files. The file path should be where csvs are on your machine or server.

importData(type = 'csv',
           filepath = "C:/NETN/R_Dev/Water/data") # update filepath to your computer

Option 2. Import data via zip file of csvs. The filepath should be the location and name of the zip file.

importData(type = 'zip',
           filepath = "C:/NETN/R_Dev/Water/data/NETN_Water_Data_Package_20240423.zip")

Option 3. Import data via data package database file on your computer

importData(type = 'dbfile',
           filepath = "C:/NETN/R_Dev/Water/data/NETN_h3Ov4_DataPackage_202331115.accdb")

Option 4. Import data via data package database DSN (Data Source Name) on your computer. Note that this is the default setting. As long as you have a named DSN called “NETNWQ_DP” that links to the data package database, and that database links to the latest NETN WQ backend database, the code below will run. See Setting up DSN tab for how to set up DSN.

importData() # easiest
importData(type = 'DSN', odbc = "NETNWQ_DP") # equivalent to line above

Step 6. (Optional) Export data package to zip

You can export all of the csvs to a zip file with the day’s date stamped on the file name. This allows you to import the tables from the database, then export the csvs as one zip file.

importData() # easiest
exportData(filepath = "./data", zip = TRUE)

Step 7. Play with the data

The functions in the waterNETN package are designed to work with the views, and are the best way to interact with the data to query by park, site, site type, year, parameter, etc. However, if you want to view the raw data, and you imported the data into the VIEWS_WQ environment, you can access them with the code below:

# See list of the views
names(VIEWS_WQ)

# View one of the views
View(VIEWS_WQ$Chemistry_Data)

# Assign a view to a data frame named chem in R. Interact with chem the way you would work with any normal data frame in R. 
chem <- VIEWS_WQ$Chemistry_Data

If you want to use the print_head() function that shows output in the markdown, run the code below. This makes the results print cleaner in the markdown report. For your purposes, you can just run: head(dataframe).

print_head <- function(df){
  knitr::kable(df[1:6,]) |> #, table.attr = "style='width:60%;'") |> 
    kableExtra::kable_classic(full_width = F, font_size = 12, 
                              bootstrap_options = c("condensed"))
}

Getting help

Getting (and improving) help

The functions in waterNETN have help documentation like any R package. To view the help, you can go to the Packages tab and click on waterNETN (see below). That will show you all the functions in the package. Clicking on individual functions will take you to the help documentation for that function.

You can also see the help of a function by running, for example:

?importData

If waterNETN isn’t loaded yet, you’d run:

?waterNETN::importData

Each function’s help includes a Description, Usage (i.e. function arguments and their defaults), Argument options/definitions, and several examples showing how the function can be used.

This is where you come in! If you notice typos or can think of better descriptions, examples, error messages, etc., please send them my way! After we’re more comfortable with R packages and get versed on GitHub, you’ll be able to make those changes directly in the package. For now, you can just send me your suggestions and I’ll make the changes.

Finally, if you ever want to peak under the hood at the function, you can view it several ways.

Keep F2 key pressed and click on the function name in R. This trick works for many but not all functions in R.
View code in the GitHub katemmiller/waterNETN repo. The functions are in the R folder.

Setting up DSN

Setting up a DSN

Go to Windows Start Menu and search ODBC. Click on ODBC Data Sources (64-bit)

Click on Add, then select Microsoft Access Driver (*.mdb, *.accdb) then click Finish in next menu.

Enter NETNWQ_DP into the Data Source Name (red arrow), click on Select (red circle). In new window, click on C:/ (orange arrow) and find the path to your database. If it’s on the Z drive, then click on the Drives window and select correct Drive. Click on the data package database (yellow arrow). When complete, click OK. Finally, add the name of the database file to the Description (blue arrow), so it’s easier to check whether you’re using the latest version.

If updating an existing DSN, follow similar process, except click on the DSN in the first window and select Configure instead of Add.

get Data functions

getClimNOAA()

Extracts gridded NOAA climate data for park centroids. After looking over Daymet, Prism and weather station data, I’ve settled on the NOAA gridded data to make comparisons because they use the same algorithm and scale to calculate historic averages and annual data. Historic averages are available for the 20th century (1901 - 2000), which I also like, and for the latest 30 year normal (1991 - 2020). The normals for NETN park centroids have been compiled and are a data source that exists when you install the waterNETN package (called NETN_clim_norms.rda). This means you don’t have to download or summarize the normals every time you want to use them, which speeds up plotting. You can download data as far back as 1951. However, I’ve also saved time by compiling all of the monthly data from Jan. 1895 - May 2024, and saved this as a dataset for the package called NETN_clim_annual.rda (see code chunk below on how to access this data set). That means you only need to download and compile data for months after May 2024. I have a webalert set for when a new month is posted, and will add that to the annual data as I have time. If I get behind, you can use the getClimNOAA() function to download the new months of data.

I’ve also updated all of the climate comparison plotting functions to work with the NOAA data. These functions will check for months not included in the NETN_clim_annual, and will download what’s available. The plotting functions use getClimNOAA() under the hood to perform this. You’re therefore unlikely to need to use the getClimNOAA() function yourself, unless you want to summarize/plot the data not using functions in waterNETN.

Currently the getClimNOAA() function only works for 1 year at a time, which should be fine, because I’ve already compiled 1895 - May 2024 data. New data availability is usually delayed by about 2 weeks after a month ends. For example, on 6/16/2024, the data for May 2024 were posted as final. If you try to download data for a month before it’s available, you’ll get an error telling you which months are not available yet. You then just need to change your months argument to include only the months that didn’t fail.

The more months you include in the function, the longer it will take, as each month and climate variable is a separate raster to download. Number of parks minimally impacts performance.

Note that I’m considering stripping most of these climate functions out of this package and creating a new climateNETN package, so any protocol can work with these data. The function names won’t change, and you won’t need to update your code beyond needing to install and load the separate package. I’ll update this website once I’m finished with that transition.

get weather data for January - May 2024 for all NETN.

This may take a minute or so to download. If you’re ultimately interested in plotting the data using one of the waterNETN package functions, don’t worry about downloading separately. Just run the function, and it will run this function under the hood for months not already compiled.

netn_jan_may <- getClimNOAA(months = 1:5)
print_head(netn_jan_may)

UNIT_CODE	long	lat	tmax	tmin	tavg	year	month
ACAD	-68.26019	44.37656	0.9296875	-3.2031250	-1.140625	2024	1
MIMA	-71.30739	42.45629	4.5390625	-1.0781250	1.726562	2024	1
SAIR	-71.00770	42.46831	3.6328125	-0.9609375	1.328125	2024	1
SARA	-73.62991	42.99831	0.8437500	-1.8437500	-0.500000	2024	1
SAGA	-72.37454	43.49917	-0.2031250	-3.2421875	-1.718750	2024	1
MORR	-74.53459	40.76522	5.2734375	1.4062500	3.343750	2024	1

Get weather data for May 2024 for ROVA

rova_may <- getClimNOAA(park = "ROVA", months = 5)
head(rova_may)

##   UnitCode
## 1     ROVA
##                                       UnitName
## 1 Roosevelt-Vanderbilt National Historic Sites
##        long
## 1 -73.92672
##        lat
## 1 41.76804
##       prcp
## 1 85.50977
##       tmax
## 1 23.55957
##       tmin
## 1 12.00977
##       tavg
## 1 17.79004
##   year
## 1 2024
##   month
## 1     5

Get weather data for all of 2023 in LNETN

lnetn_2023 <- getClimNOAA(park = "LNETN", year = 2023, months = 1:12)
print_head(lnetn_2023)

UNIT_CODE	long	lat	prcp	tmax	tmin	tavg	year	month
MIMA	-71.30739	42.45629	6.000000	13.53906	5.8906250	9.710938	2023	1
SAIR	-71.00770	42.46831	7.648438	13.09375	6.5312500	9.812500	2023	1
SARA	-73.62991	42.99831	1.523438	9.43750	1.9531250	5.687500	2023	1
SAGA	-72.37454	43.49917	1.296875	8.34375	0.7265625	4.539062	2023	1
MORR	-74.53459	40.76522	7.398438	10.82031	5.4062500	8.117188	2023	1
WEFA	-73.45453	41.25886	8.179688	12.57812	5.9609375	9.273438	2023	1

View dataset containing 1895 - 2024 NOAA climate data

data("NETN_clim_annual")
print_head(NETN_clim_annual)

UnitCode	UnitName	long	lat	prcp	tmax	tmin	tavg	year	month
ACAD	Acadia National Park	-68.26019	44.37656	150.50000	-0.5898438	-12.139648	-6.370117	1895	1
BOHA	Boston Harbor Islands National Recreation Area	-71.02159	42.27097	96.19043	0.7099609	-8.129883	-3.709961	1895	1
ROVA	Roosevelt-Vanderbilt National Historic Sites	-73.92672	41.76804	92.76953	-1.7099609	-10.780273	-6.250000	1895	1
MABI	Marsh-Billings-Rockefeller National Historical Park	-72.53824	43.63382	55.83984	-2.6503906	-16.750000	-9.700195	1895	1
MIMA	Minute Man National Historical Park	-71.29604	42.45340	92.24023	0.1298828	-10.809570	-5.339844	1895	1
MORR	Morristown National Historical Park	-74.53749	40.76780	118.08984	0.2197266	-9.500000	-4.639648	1895	1

getClimDaymet()

NOTE THAT THIS IS NO LONGER THE PRIMARY FUNCTION TO USE FOR WORKING WITH CLIMATE DATA. . I’m leaving it in the package just in case there’s a reason to use it in the future. If I create a new climateNETN package, I’ll remove this function from the waterNETN package.

Daymet data are 1km gridded daily climate data, similar to PRISM data, but with more daily data associated with each location and more efficient downloading process (extract data from lat/long points, rather than download nationwide raster files and extract manually). Daymet data go back to 1980. See https://daymet.ornl.gov/ for more information. Another benefit of Daymet over PRISM is the data can be used to calculate drought indices, like Standardized Precipitation-Evapotranspiration Index (SPEI). Note: I’ll eventually add SPEI metrics to the output of this function.

Daymet data are typically available for a calendar year a a few months into the next year.

In waterNETN, you can download Daymet data by park, site, site type, and spanning years 1980 to 2023. See examples below. The larger the number of sites and years, the longer it can take.

The following code downloads Daymet data for East Primrose Brook in MORR for the entire monitoring period.

prim <- getClimDaymet(site = "MORRSB", years = 2006:2023)
print_head(prim) # table of first 6 rows- this is a custom function I wrote to save on coding

The following code downloads Daymet data for all sites in MABI for 2023. Results are assigned to R session environment and not exported to data folder in the working directory.

# Daymet data for the Pogue and Pogue stream for 2023, but not exporting to disk.
mabi_dm <- getClimDaymet(park = "MABI", years = 2023)
print_head(mabi_dm)

SiteCode	dm_tile	SiteLatitude	SiteLongitude	altitude	Date	year	yday	dm_dayl_s	dm_prcp_mmday	dm_srad_Wm2	dm_tmax_degc	dm_tmin_degc	dm_vp_Pa
MABISA	11934	43.63493	-72.52937	302	2023-01-02	2023	1	31737.91	0.00	166.06	7.88	0.32	624.97
MABISA	11934	43.63493	-72.52937	302	2023-01-03	2023	2	31785.76	0.00	169.12	6.08	-1.57	544.55
MABISA	11934	43.63493	-72.52937	302	2023-01-04	2023	3	31837.50	5.54	105.08	2.20	-3.49	472.10
MABISA	11934	43.63493	-72.52937	302	2023-01-05	2023	4	31893.11	11.65	110.11	4.30	-1.77	536.58
MABISA	11934	43.63493	-72.52937	302	2023-01-06	2023	5	31952.54	2.13	71.89	4.10	0.43	630.18
MABISA	11934	43.63493	-72.52937	302	2023-01-07	2023	6	32015.75	4.29	45.28	2.27	0.09	614.84

getClimWStat()

For each NETN monitoring site, I used the rnoaa package to find all weather stations within 20km radius of each site’s lat/long coordinates. I then checked the period of record for each of these sites. I selected the closest weather station with the longest period of record that included 2023 and 2024. This assessment was done at the park level, so the data for each site within a park is the same. The API used to download weather station data is data.rcc-acis.org. Unfortunately only temperature data are available for ACAD’s McFarland Hill weather station (AKA MARS). ACAD’s hourly precipitation data were downloaded separately from the NADP website http://nadp2.slh.wisc.edu/siteOps/ppt/. These data only go back to 2008 and were aggregated to daily total precipitation. This function is slower for ACAD as a result. Though note that data are downloaded at the park level, then joined with site-specific data. In other words, data are downloaded for a park x combination of years once, then joined to site-level data. Weather stations used and their distance from each site are included in the output.

The following code downloads weather station data for all sites in ROVA the past 10 years and writes to a csv.

rova_ws <- getClimWStat(park = "ROVA", years = 2013:2023, export = T, 
                        filepath = "./data")
print_head(rova_ws)

SiteCode	SiteName	UnitCode	Date	year	month	doy	ws_id	ws_lat	ws_long	ws_dist_km	ws_tmaxc	ws_tminc	ws_pcpmm
ROVASA	Lower FDR Brook	ROVA	2013-01-01	2013	1	1	USW00064756	41.7858	-73.7422	16.564	1.667	-10.000	0.000
ROVASA	Lower FDR Brook	ROVA	2013-01-02	2013	1	2	USW00064756	41.7858	-73.7422	16.564	-1.667	-17.222	0.000
ROVASA	Lower FDR Brook	ROVA	2013-01-03	2013	1	3	USW00064756	41.7858	-73.7422	16.564	-5.000	-18.333	0.000
ROVASA	Lower FDR Brook	ROVA	2013-01-04	2013	1	4	USW00064756	41.7858	-73.7422	16.564	2.778	-6.667	0.000
ROVASA	Lower FDR Brook	ROVA	2013-01-05	2013	1	5	USW00064756	41.7858	-73.7422	16.564	3.889	-6.667	0.000
ROVASA	Lower FDR Brook	ROVA	2013-01-06	2013	1	6	USW00064756	41.7858	-73.7422	16.564	5.000	-5.556	0.254

The following code downloads weather station data for all lakes in ACAD in 2023, but doesn’t export to disk.

acad_ws <- getClimWStat(park = "ACAD", site_type = "lake", years = 2023)
print_head(acad_ws)

SiteCode	SiteName	UnitCode	Date	year	month	doy	ws_id	ws_lat	ws_long	ws_dist_km	ws_tmaxc	ws_tminc	ws_pcpmm
ACANTB	Aunt Bettys Pond	ACAD	2023-01-01	2023	1	1	USR0000MMCF	44.3769	-68.2608	1.319	11.667	1.111	0.000
ACANTB	Aunt Bettys Pond	ACAD	2023-01-02	2023	1	2	USR0000MMCF	44.3769	-68.2608	1.319	6.667	0.000	0.000
ACANTB	Aunt Bettys Pond	ACAD	2023-01-03	2023	1	3	USR0000MMCF	44.3769	-68.2608	1.319	5.556	1.111	3.556
ACANTB	Aunt Bettys Pond	ACAD	2023-01-04	2023	1	4	USR0000MMCF	44.3769	-68.2608	1.319	5.556	-0.556	4.064
ACANTB	Aunt Bettys Pond	ACAD	2023-01-05	2023	1	5	USR0000MMCF	44.3769	-68.2608	1.319	-0.556	-2.778	5.588
ACANTB	Aunt Bettys Pond	ACAD	2023-01-06	2023	1	6	USR0000MMCF	44.3769	-68.2608	1.319	-1.111	-2.778	3.556

getClimDrought()

This function downloads weekly drought index at the county level based on the U.S. Drought Monitor for each specified park. If downloading for multiple parks and multiple weeks/years, function may be slow. Returned data frame includes percent of county area in 5 levels of drought, with D0 = Abnormally Dry, D1 = Moderate Drought, D2 = Severe Drought, D3 = Extreme Drought, and D4 = Exceptional Drought. Also returned is the Drought Severity and Coverage Index, which ranges from 0 to 500, and is a weighted sum of area within each of the drought categories. A score of 500 indicates the entire area is in exceptional drought (D4). Where multiple counties occur in a park, data can be faceted by county if dom_county = FALSE. To only plot predominant county, specify dom_county = TRUE (default). Otherwise, only park-level data are returned, as the results won’t vary by site within a park. Note that Drought Monitor reports data such that, when drought a level moves into a higher category (e.g. D1 = 10), the previous level (e.g. D0) is listed as 100. The columns with pct (e.g. D0pct) correct for that, so that D0 + … + D05 + None = 100.For more info on the data, see <US Drought Monitor

Get drought info for MABI and SAGA sites for first week of May

mabisaga <- getClimDrought(park = c("MABI", "SAGA"), week_start = "05/01/2024")
print_head(mabisaga)

	MapDate	FIPS	County	State	None	D0	D1	D2	D3	D4	ValidStart	ValidEnd	StatisticFormatID	UnitCode	DSCI	D0pct	D1pct	D2pct	D3pct	D4pct
1	20240507	50027	Windsor County	VT	100	0	0	0	0	0	2024-05-07	2024-05-13	1	MABI	0	0	0	0	0	0
2	20240430	50027	Windsor County	VT	100	0	0	0	0	0	2024-04-30	2024-05-06	1	MABI	0	0	0	0	0	0
3	20240507	33019	Sullivan County	NH	100	0	0	0	0	0	2024-05-07	2024-05-13	1	SAGA	0	0	0	0	0	0
4	20240430	33019	Sullivan County	NH	100	0	0	0	0	0	2024-04-30	2024-05-06	1	SAGA	0	0	0	0	0	0
NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.1	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA

Get drought info for both counties in MORR in 2023

morr <- getClimDrought(park = "MORR", years = 2023, dom_county = TRUE)
print_head(morr)

MapDate	FIPS	County	State	None	D0	ValidStart	ValidEnd	StatisticFormatID	UnitCode	DSCI	D0pct
20231226	34027	Morris County	NJ	100.00	0.00	2023-12-26	2024-01-01	1	MORR	0	0.00
20231219	34027	Morris County	NJ	100.00	0.00	2023-12-19	2023-12-25	1	MORR	0	0.00
20231212	34027	Morris County	NJ	100.00	0.00	2023-12-12	2023-12-18	1	MORR	0	0.00
20231205	34027	Morris County	NJ	64.42	35.58	2023-12-05	2023-12-11	1	MORR	36	35.58
20231128	34027	Morris County	NJ	64.42	35.58	2023-11-28	2023-12-04	1	MORR	36	35.58
20231121	34027	Morris County	NJ	9.71	90.29	2023-11-21	2023-11-27	1	MORR	90	90.29

getSites()

Query site-level data. This function combines columns in common between the Sites_Stream and Sites_Lake views. This is a good building block for other functions, but may be less helpful on its own. Though, one helpful use of this function is to get the site codes for a given park or site.

Get site info for ACAD lakes

ACAD_lakes <- getSites(park = "ACAD", site_type = "lake")

Get site info for the Pogue

pogue <- getSites(site = "MABIPA")

List site codes and their full site name for SARA

sara_sites <- getSites(park = "SARA") |> select(SiteCode, SiteName)
print_head(sara_sites)

	SiteCode	SiteName
1	SARASA	Kroma Kill
2	SARASC	Upper Mill Creek
3	SARASD	Mill Creek Confluence
NA	NA	NA
NA.1	NA	NA
NA.2	NA	NA

getSitesLake()

Query site-level data for lakes only. This returns lake-only columns in the Sites_Lake views

Get site info for ACAD lakes

ACAD_str <- getSitesLake(park = "ACAD")
print_head(ACAD_str)

SiteCode	SiteName	UnitCode	SubUnitCode	HUC12	ContribWshedArea_km2	Inlet	Outlet	MaxDepth_m	TrophicType	LegislativeClass	FishMgtType	SiteDescription	Notes
ACANTB	Aunt Bettys Pond	ACAD	NA	10500021502	2.77	unnamed	RICHARDSON BRK	2.1	EUTROPHIC	GPA	COLDWATER	Shallow, marshy pond with an area of 31.5 acres. Inlet: Aunt Betty Inlet. Outlet: Richardson Brook.	Pond is full of emergent vegetation by mid-summer. Deep hole is in NE corner, near outlet. Surrounded by ANP. Last surveyed by IFW in 1982. Closed to ice fishing.
ACBOWL	The Bowl	ACAD	NA	10500021501	0.21	NA	OTTER CREEK	9.0	MESOTROPHIC	GPA	NONE	High elevation, mesotrophic pond located between Cahamplaine Mtn. and The Beehive.	Located off the Beehive Trail. Fish survey conducted in 1995 BY J.R. Burgess. Not monitored by NETN- sampled annually for acid depositon effects by ACAD staff. Once served as water supply for Satterlee estate.
ACBRBK	Bear Brook Pond	ACAD	NA	10500021501	NA	BEAR BROOK	BEAR BROOK	3.9	MESOTROPHIC	GPA	NONE	Mesotrophic, unstratified pond with an area of 7.5 acres. Its inlet and outlet is Bear Brook.	Located on Park Loop Road near Jackson Lab. Site of ANP purple loosestrife monitoring/control plots. Active beaver in area. Fish survey in 1995 by J.R. Burgess. Also known as Beaver Dam Pond.
ACBUBL	Bubble Pond	ACAD	NA	10500021501	1.77	unnamed	BUBBLE BROOK	11.9	MESOTROPHIC	GPA	COLDWATER	Oligotrophic, stratified lake with an area of 33 acres. Part of the public water supply. Inlet: unnamed stream on southern end. Outlet: Bubble Brook.	Surrounded by ANP. Connects with Bar Harbor water supply. Last surveyed by IFW in 1942. Closed to ice fishing.
ACEAGL	Eagle Lake	ACAD	NA	10500021501	5.60	BUBBLE BROOK	DUCK BROOK	33.5	OLIGOTROPHIC	GPA	COLDWATER	Oligotrophic, stratified lake with an area of 465.8 acres. Part of the public water supply. Inlet: Bubble Brook. Outlet: Duck Brook.	LLS spawning area built in outlet in 55. 10 hp outboard limit. Last surveyed by IFW in 1990. Surrounded by ANP.
ACECHO	Echo Lake	ACAD	NA	10500021502	5.10	LURVEY SPRG BRK	DENNING BROOK	20.1	OLIGOTROPHIC	GPA	COLDWATER	Oligotrophic, stratified lake with an area of 236.8 acres. Inlet: Lurvey Spring Brook. Outlet: Denning Brook.	Partial boundary w/ANP. Last surveyed by IFW in 1988. Outboard hp restrictions. Reclaimed in 1956. Guarded ANP swim beach on south end.

Get all site info for the Pogue Pond

pogue <- getSitesLake(site = "MABIPA", output = 'verbose')
print_head(pogue)

	GroupCode	GroupName	UnitCode	UnitName	SubUnitCode	SubUnitName	SiteCode	SiteName	SiteLatitude	SiteLongitude	Datum	XYAccuracy	HUC12	ContribWshedArea_km2	Inlet	Outlet	MaxDepth_m	TrophicType	LegislativeClass	FishMgtType	SiteDescription	Notes	IsPointCUI
1	NETN	Northeast Temperate Network	MABI	Marsh-Billings-Rockefeller National Historical Park	NA	NA	MABIPA	The Pogue	43.63336	-72.54263	NAD83	NA	10801060206	0.57	NA	Pogue Brk	3.4	MESOTROPHIC	A(1)	COLDWATER	The Pogue is a 14-acre pond at the headwaters of the Pogue Brook. Naturally a spring-fed boggy area, it was created in the 1880s with an earthen dam. Swimming, fishing, or wading is prohibited.		FALSE
NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.1	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.2	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.3	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.4	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA

getSitesStream()

Query site-level data for streams only. This returns stream-only columns in the Sites_Stream views

Get site info for ACAD streams

ACAD_str <- getSitesStream(park = "ACAD")
print_head(ACAD_str)

SiteCode	UnitCode	SubUnitCode	USGS_StaNumb	ContribWshedArea_km2	SiteDescription	Notes	LegislativeClass	FisheryType
ACABIN	ACAD	NA	1022869	2.4802373	The inlet flows north from Gilmore Meadow into Aunt Betty Pond.	Beavers often present both US & DS of carriage road.	AA	Warmwater
ACBRKB	ACAD	NA	1022825	3.7293399	The stream flows north from the Breakneck Ponds and into the ocean (Hulls Cove).	Exact sampling point can vary depending on streamflow/water level. Hits baseflow in dry years.	AA	Coldwater
ACBRWN	ACAD	NA	1022866	1.1840292	The stream flows northwest from Parkman Mountain into the ocean (Somes Sound). Measure east (US) of culvert under Sargent Drive.	Very close to municipal sand/salt storage facility.	AA	Coldwater
ACCADS	ACAD	NA	1022835	0.6223256	On right bank approximately 500 feet upstream of the crossing of Canon Brook Trail (500 ft upstream from confluence with Otter Creek and 0.5 miles southeast of Cadillac summit)	High-elevation stream flows south on the east face of Cadillac mountain into Canon Brook- Otter Creek. Former site of USGS stream gage 01022835 (1999- 2006).	AA	Coldwater
ACDKLI	ACAD	NA	10228755	0.4818996	The stream flows north from Duck Pond into Long Pond (MDI).	Can be near or at baseflow in dry years and difficult to measure streamflow.	AA	Coldwater
ACDUCK	ACAD	NA	1022827	9.6264998	Sample site is located north of Rte. 233 crossing, approx. 200 ft. DS of culvert.	Duck Brook flows north from Eagle Lake into the ocean (Hulls Cove/Frenchman Bay). Established in 2009, replaces ACEGLO.	AA	Coldwater

Get all site info for the Pogue Stream

poguestr <- getSitesStream(site = "MABISA", output = 'verbose')
print_head(poguestr)

	GroupCode	GroupName	UnitCode	UnitName	SubUnitCode	SubUnitName	SiteCode	SiteName	SiteLatitude	SiteLongitude	Datum	XYAccuracy	USGS_StaNumb	ContribWshedArea_km2	SiteDescription	Notes	LegislativeClass	FisheryType	IsPointCUI
1	NETN	Northeast Temperate Network	MABI	Marsh-Billings-Rockefeller National Historical Park	NA	NA	MABISA	Pogue Brook	43.63493	-72.52937	NAD83	NA	NA	0.9945554	Only stream within the park, the outlet of The Pogue. One stream-sampling site was selected to represent Pogue Brook. This site was chosen because it is the most downstream location within park boundaries.		A(1)	Coldwater	FALSE
NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.1	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.2	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.3	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA
NA.4	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA	NA

getEvents()

Query water data by park, site, site_type, year and month.

Get event info for Witch Hole Pond in ACAD in 2023

ACAD_lakes <- getEvents(site = "ACWHOL", years = 2023)
print_head(ACAD_lakes)

UnitCode	SiteCode	SiteName	SiteType	SiteLatitude	SiteLongitude	EventDate	EventCode	year	month	doy
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-05-03	4649	2023	5	122
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-05-17	4652	2023	5	136
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-06-26	4698	2023	6	176
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-07-26	4741	2023	7	206
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-08-21	4782	2023	8	232
ACAD	ACWHOL	Witch Hole Pond	Lake	44.40001	-68.24295	2023-09-18	4818	2023	9	260

Get event info for all years in the Pogue

pogue <- getEvents(site = "MABIPA", years = 2006:2023)
print_head(pogue)

UnitCode	SiteCode	SiteName	SiteType	SiteLatitude	SiteLongitude	EventDate	EventCode	year	month	doy
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2006-06-13	1330	2006	6	163
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2006-07-20	671	2006	7	200
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2006-08-23	1331	2006	8	234
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2006-09-13	1054	2006	9	255
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2006-10-11	1055	2006	10	283
MABI	MABIPA	The Pogue	Lake	43.63336	-72.54263	2007-05-08	1057	2007	5	127

Get date of sampling each site in lower NETN in May, 2023

lnetn_evs <- getEvents(park = "LNETN", years = 2023, months = 5) |> 
  select(SiteCode, SiteName, EventDate)
print_head(lnetn_evs) #just printing first 6 rows

SiteCode	SiteName	EventDate
MABIPA	The Pogue	2023-05-16
MABISA	Pogue Brook	2023-05-16
MIMASA	Mill Brook	2023-05-02
MIMASB	Elm Brook	2023-05-02
MIMASC	Concord River	2023-05-02
MORRSB	Primrose Brook Confluence	2023-05-10

getChemistry()

This function allows you to query the Chemistry_Data view by park, site, site type, year, month, parameter, and sample type. The returned data frame is long (i.e. stacked) to facilitate data summary and plotting. Note that sample depth is not a filter in getChemistry() as it is with the Sonde In Situ data.

Get N parameters for all sites and non-QAQC events in MIMA.

Note that QC_type = "ENV" is the default for this function, which returns only non-QAQC events. Note also the use of named objects for the arguments. This allows you to set them at the top of a script, rather than having to type them out repeatedly. You can then change them in 1 place (i.e., update the year to 2024) and rerun the code.

n_params <- c("NH3", "NH3_mgL", "NO2", "NO2_mgL", "NO2+NO3", 
              "NO2+NO3_mgL", "NO3", "NO3_ueqL", "TN", "TN_mgL")
period <- 2006:2023
mima_n <- getChemistry(park = "MIMA", years = period, parameter = n_params)
print_head(mima_n) # top 6 rows

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	datetime	QCtype	SampleType	SampleTime	SampleDepth_m	Parameter	Value	ValueFlag	ValueUnit	ValueDetectionCondition	censored	LabCode
MIMASA	Mill Brook	MIMA	NA	2006-08-21	766	2006	8	232	2006-08-21 12:45:00	ENV	G	12:45:00	0.2	NO2+NO3_mgL	0.169	NA	mg/L	Detected And Quantified	FALSE	NWQL
MIMASA	Mill Brook	MIMA	NA	2006-08-21	766	2006	8	232	2006-08-21 12:45:00	ENV	G	12:45:00	0.2	NO2_mgL	0.002	NA	mg/L	Detected And Quantified	FALSE	NWQL
MIMASA	Mill Brook	MIMA	NA	2006-08-21	766	2006	8	232	2006-08-21 12:45:00	ENV	G	12:45:00	0.2	NH3_mgL	0.030	NA	mg/L	Detected And Quantified	FALSE	NWQL
MIMASA	Mill Brook	MIMA	NA	2006-08-21	766	2006	8	232	2006-08-21 12:45:00	ENV	G	12:45:00	0.2	TN_mgL	0.560	NA	mg/L	Detected And Quantified	FALSE	NWQL
MIMASA	Mill Brook	MIMA	NA	2007-05-07	770	2007	5	126	2007-05-07 12:31:00	ENV	G	12:31:00	0.2	NO2+NO3_mgL	0.338	NA	mg/L	Detected And Quantified	FALSE	NWQL
MIMASA	Mill Brook	MIMA	NA	2007-05-07	770	2007	5	126	2007-05-07 12:31:00	ENV	G	12:31:00	0.2	TN_mgL	0.840	NA	mg/L	Detected And Quantified	FALSE	NWQL

Get results for all surface samples, including QC samples, in Browns Brook in 2023 in ACAD

brbk <- getChemistry(site = "ACBRWN", QC_type = "all", param = "all", 
                     years = 2023)
print_head(brbk)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	datetime	QCtype	SampleType	SampleTime	SampleDepth_m	Parameter	Value	ValueFlag	ValueUnit	ValueDetectionCondition	censored	LabCode
ACBRWN	Browns Brook	ACAD	NA	2023-05-08	4669	2023	5	127	2023-05-08 15:32:47	ENV	G	15:32:47	0.1	TP_ugL	1.500	NA	ug/L	Detected And Quantified	FALSE	SWRL
ACBRWN	Browns Brook	ACAD	NA	2023-05-08	4669	2023	5	127	2023-05-08 15:32:47	ENV	G	15:32:47	0.1	DOC_mgL	0.900	NA	mg/L	Detected And Quantified	FALSE	SWRL
ACBRWN	Browns Brook	ACAD	NA	2023-05-08	4669	2023	5	127	2023-05-08 15:32:47	ENV	G	15:32:47	0.1	pH_Lab	6.660	NA	pH std units	Detected And Quantified	FALSE	SWRL
ACBRWN	Browns Brook	ACAD	NA	2023-05-08	4669	2023	5	127	2023-05-08 15:32:47	ENV	G	15:32:47	0.1	TN_mgL	0.038	NA	mg/L	Detected And Quantified	FALSE	SWRL
ACBRWN	Browns Brook	ACAD	NA	2023-05-08	4669	2023	5	127	2023-05-08 15:32:47	ENV	G	15:32:47	0.1	ANC_ueqL	57.600	NA	ueq/L	Detected And Quantified	FALSE	SWRL
ACBRWN	Browns Brook	ACAD	NA	2023-06-08	4710	2023	6	158	2023-06-08 14:03:51	ENV	G	14:03:51	0.1	pH_Lab	6.540	NA	pH std units	Detected And Quantified	FALSE	SWRL

Get ANC for LNETN parks for all years during month of May. Note that by not specifying years, all years by default will be included in the output.

anc <- getChemistry(park = "LNETN", param = "ANC", months = 5)
print_head(anc)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	datetime	QCtype	SampleType	SampleTime	SampleDepth_m	Parameter	Value	ValueFlag	ValueUnit	ValueDetectionCondition	censored	LabCode
MABISA	Pogue Brook	MABI	NA	2007-05-08	676	2007	5	127	2007-05-08 13:25:00	ENV	G	13:25:00	0.2	ANC_ueqL	1416.2	NA	ueq/L	Detected And Quantified	FALSE	SECRL
MABISA	Pogue Brook	MABI	NA	2008-05-20	1107	2008	5	140	NA	ENV	G	NA	0.2	ANC_ueqL	1435.0	NA	ueq/L	Detected And Quantified	FALSE	SECRL
MABISA	Pogue Brook	MABI	NA	2009-05-07	1224	2009	5	126	NA	ENV	G	NA	0.2	ANC_ueqL	1370.0	NA	ueq/L	Detected And Quantified	FALSE	SECRL
MABISA	Pogue Brook	MABI	NA	2010-05-25	236	2010	5	144	NA	ENV	G	NA	0.2	ANC_ueqL	1770.0	NA	ueq/L	Detected And Quantified	FALSE	SECRL
MABISA	Pogue Brook	MABI	NA	2011-05-24	1747	2011	5	143	2011-05-24 14:30:00	ENV	G	14:30:00	0.2	ANC_ueqL	1390.0	NA	ueq/L	Detected And Quantified	FALSE	SECRL
MABISA	Pogue Brook	MABI	NA	2012-05-14	1977	2012	5	134	2012-05-14 10:20:00	ENV	G	10:20:00	0.2	ANC_ueqL	1570.0	NA	ueq/L	Detected And Quantified	FALSE	SECRL

Get for ACAD lakes from 2021:2023 from all sample depths.

surf <- getChemistry(park = "ACAD", site_type = "lake", years = 2021:2023)
print_head(surf)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	datetime	QCtype	SampleType	SampleTime	SampleDepth_m	Parameter	Value	ValueFlag	ValueUnit	ValueDetectionCondition	censored	LabCode
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-05-24	4171	2021	5	143	2021-05-24 10:26:22	ENV	G	10:26:22	0.2	TP_ugL	11.000	NA	ug/L	Detected And Quantified	FALSE	SWRL
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-05-24	4171	2021	5	143	2021-05-24 10:26:22	ENV	G	10:26:22	0.2	TN_mgL	0.372	NA	mg/L	Detected And Quantified	FALSE	SWRL
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-05-24	4171	2021	5	143	2021-05-24 10:26:22	ENV	G	10:26:22	0.2	DOC_mgL	5.100	NA	mg/L	Detected And Quantified	FALSE	SWRL
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-05-24	4171	2021	5	143	2021-05-24 10:26:22	ENV	G	10:26:22	0.2	ChlA_ugL	2.300	NA	ug/L	Detected And Quantified	FALSE	SWRL
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-06-24	4197	2021	6	174	2021-06-24 11:07:28	ENV	G	11:07:28	0.2	TP_ugL	13.100	NA	ug/L	Detected And Quantified	FALSE	SWRL
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-06-24	4197	2021	6	174	2021-06-24 11:07:28	ENV	G	11:07:28	0.2	TN_mgL	0.462	NA	mg/L	Detected And Quantified	FALSE	SWRL

Get censored and non-censored NH3 data for LNETN all years. Note the Flag column indicates the measurement is censored by reporting the detection limit used.

cens <- getChemistry(park = "LNETN", param = "NH3", include_censored = TRUE)
print_head(cens)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	datetime	QCtype	SampleType	SampleTime	SampleDepth_m	Parameter	Value	ValueFlag	ValueUnit	ValueDetectionCondition	censored	LabCode
MABIPA	The Pogue	MABI	NA	2006-06-13	1330	2006	6	163	2006-06-13 15:30:00	ENV	G	15:30:00	0.2	NH3_mgL	0.006	E, <MRL 0.01	mg/L	Detected And Quantified	FALSE	NWQL
MABIPA	The Pogue	MABI	NA	2006-06-13	1330	2006	6	163	2006-06-13 15:30:00	ENV	G	15:30:00	0.2	NH3_mgL	0.010	E, <MRL 0.01	mg/L	NA	TRUE	NWQL
MABIPA	The Pogue	MABI	NA	2006-08-23	1331	2006	8	234	2006-08-23 10:45:00	ENV	G	10:45:00	0.2	NH3_mgL	0.068	NA	mg/L	Detected And Quantified	FALSE	NWQL
MABIPA	The Pogue	MABI	NA	2007-06-07	1058	2007	6	157	2007-06-07 14:30:00	ENV	G	14:30:00	0.2	NH3_mgL	0.053	NA	mg/L	Detected And Quantified	FALSE	NWQL
MABIPA	The Pogue	MABI	NA	2007-08-14	1060	2007	8	225	2007-08-14 14:20:00	ENV	G	14:20:00	0.2	NH3_mgL	0.099	NA	mg/L	Detected And Quantified	FALSE	NWQL
MABIPA	The Pogue	MABI	NA	2008-06-17	1101	2008	6	168	NA	ENV	G	NA	0.2	NH3_mgL	0.080	<MRL 0.08	mg/L	NA	TRUE	SECRL

getDischarge()

This function allows you to query the Discharge_Data view by park, site, site type, year, month, measurement method (e.g., flowtracker or pygmy), and measurement rating (e.g., E, G, F, P). Note that you can also return all columns or a reduced set of columns with the output argument. Default is output = ‘short’. This function is set up to work with site_type = ‘stream’, so you don’t have to specify that.

Get discharge data for all LNETN sites and years measured with flowtracker

ft <- getDischarge(park = "LNETN", method = "Pygmy")
print_head(ft)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	ReachType	FlowStatus	DischargeMethod	TotalArea_sqft	AvgVel_fs	VelocityFlag	Discharge_cfs	DischargeFlag	MeasurementRating	Comments
MABISA	Pogue Brook	MABI	NA	2006-05-16	669	2006	5	135	NA	NA	Pygmy	1.98	1.06	NA	2.50	NA	NA	NA
MABISA	Pogue Brook	MABI	NA	2006-06-13	670	2006	6	163	Riffle	Stable	Pygmy	1.21	0.65	NA	0.95	NA	NA	post 99 DO level was low, recalibrated and retook 0 and post 99
MABISA	Pogue Brook	MABI	NA	2006-07-20	779	2006	7	200	Chute	Stable	Pygmy	0.68	0.55	NA	0.35	NA	NA	NA
MABISA	Pogue Brook	MABI	NA	2006-08-23	672	2006	8	234	Run	Stable	Pygmy	0.31	0.45	NA	0.13	NA	NA	no satellite signal
MABISA	Pogue Brook	MABI	NA	2006-09-13	673	2006	9	255	Chute	Stable	Pygmy	0.29	0.35	NA	0.11	NA	NA	NA
MABISA	Pogue Brook	MABI	NA	2006-10-11	674	2006	10	283	Run	Stable	Pygmy	0.48	0.23	NA	0.11	NA	NA	NA

Get discharge data for all ACAD streams with measurement rating of Excellent or Good

acad <- getDischarge(park = "ACAD", rating = c("E", "G"))
str(acad) # structure of returned data

## 'data.frame':    156 obs. of  19 variables:
##  $ SiteCode         : chr  "ACABIN" "ACABIN" "ACABIN" "ACABIN" ...
##  $ SiteName         : chr  "Aunt Betty Pond Inlet" "Aunt Betty Pond Inlet" "Aunt Betty Pond Inlet" "Aunt Betty Pond Inlet" ...
##  $ UnitCode         : chr  "ACAD" "ACAD" "ACAD" "ACAD" ...
##  $ SubUnitCode      : chr  NA NA NA NA ...
##  $ EventDate        : POSIXct, format:  ...
##  $ EventCode        : int  3541 3584 3775 4076 4130 4400 4485 4573 3540 3642 ...
##  $ year             : num  2018 2018 2018 2020 2020 ...
##  $ month            : num  5 6 10 7 9 5 7 9 5 7 ...
##  $ doy              : num  126 157 282 196 253 138 186 243 121 185 ...
##  $ ReachType        : chr  "Run" "Run" "Run" "Run" ...
##  $ FlowStatus       : chr  "Stable" "Stable" "Stable" "Stable" ...
##  $ DischargeMethod  : chr  "Flowtracker" "Flowtracker" "Flowtracker" "Flowtracker" ...
##  $ TotalArea_sqft   : num  20.21 13.03 18.2 10.19 6.51 ...
##  $ AvgVel_fs        : num  0.06 0.02 0.03 0.02 0 0.06 0.26 0.19 0.41 0.17 ...
##  $ VelocityFlag     : chr  NA NA NA NA ...
##  $ Discharge_cfs    : num  1.21 0.21 0.56 0.18 0.01 0.76 0.46 1.23 1.81 0.25 ...
##  $ DischargeFlag    : chr  NA NA NA NA ...
##  $ MeasurementRating: chr  "G" "G" "G" "G" ...
##  $ Comments         : chr  "Completed 2pt readings, water level higher then 1.5\nKeep file with no operator initials " "QM was taken between two beaver dams.  Flow level is bank full in this section because of the placement of the "| __truncated__ "Pool in some areas along the cross section " "Beaver dam downstream of site. Photos attached" ...

getLightPen()

This function allows you to query the Light_Penetration_Data view by park, site, site type, year, and month. Again output = ‘short’ is the default. If specify output = ‘verbose’, you’ll get all possible columns. This is only set up to work with site_type = ‘lake’, so you don’t have to specify that.

Get light penetration for Weir Pond from 2013 - 2023 all months.

period <- 2013:2023
wefa <- getLightPen(park = "WEFA", years = period)
print_head(wefa) # top 6 rows

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	MeasurementTime	MeasurementDepth_m	LightDeck	LightUW	PenetrationRatio
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:24:59	0.10	567.56	561.930	0.99010
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:25:28	0.25	564.55	454.740	0.80551
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:26:08	0.50	572.55	353.580	0.61777
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:26:59	0.75	698.19	306.090	0.43824
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:27:48	1.00	1576.10	466.950	0.29163
WEFAPA	Weir Pond	WEFA	NA	2013-05-14	2227	2013	5	133	13:28:29	1.25	702.95	141.160	0.20073

Get light penetration for ACAD lakes in 2023 in July and August.

acad <- getLightPen(park = "ACAD", years = 2023, months = 7:8)
print_head(acad) # top 6 rows

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	MeasurementTime	MeasurementDepth_m	LightDeck	LightUW	PenetrationRatio
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:37:03	0.1	654.66	635.000	0.969970
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:37:59	0.5	694.73	434.830	0.625900
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:38:45	1.0	760.40	410.100	0.539320
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:52:46	10.0	1080.70	26.797	0.024795
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:54:55	11.0	1117.10	21.420	0.019174
ACJORD	Jordan Pond	ACAD	NA	2023-07-19	4736	2023	7	199	11:58:58	12.0	1004.10	11.948	0.011899

getSecchi()

This function allows you to query Secchi data by park, site, site type, years, months and whether you want to return all observations, the first or the second. Default returns all observations with data. Note that this function is designed to only work with lakes, and doesn’t require site_type = “lake” to be specified.

Get Secchi depth the Pogue from 2021-2023, first observer only

mabi <- getSecchi(site = "MABIPA", years = 2021:2023, observer_type = "first")
print_head(mabi)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	Parameter	Value	SD_HitBottom	Observer
MABIPA	The Pogue	MABI	NA	2021-05-24	4178	2021	5	143	SecchiDepth_m	2.80	TRUE	1
MABIPA	The Pogue	MABI	NA	2021-06-22	4192	2021	6	172	SecchiDepth_m	3.00	TRUE	1
MABIPA	The Pogue	MABI	NA	2021-07-20	4239	2021	7	200	SecchiDepth_m	2.90	TRUE	1
MABIPA	The Pogue	MABI	NA	2021-08-25	4281	2021	8	236	SecchiDepth_m	3.01	TRUE	1
MABIPA	The Pogue	MABI	NA	2021-09-21	4315	2021	9	263	SecchiDepth_m	3.05	TRUE	1
MABIPA	The Pogue	MABI	NA	2021-10-22	4349	2021	10	294	SecchiDepth_m	2.85	TRUE	1

Get Secchi for all ACAD lakes sampled in July for all observers

ACAD_lake <- getSecchi(park = 'ACAD', months = 7)
print_head(ACAD_lake)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	Parameter	Value	SD_HitBottom	Observer
ACANTB	Aunt Bettys Pond	ACAD	NA	2006-07-13	38	2006	7	193	SecchiDepth_m	1.95	TRUE	1
ACANTB	Aunt Bettys Pond	ACAD	NA	2009-07-23	1136	2009	7	203	SecchiDepth_m	2.34	TRUE	1
ACANTB	Aunt Bettys Pond	ACAD	NA	2012-07-24	1389	2012	7	205	SecchiDepth_m	2.44	FALSE	1
ACANTB	Aunt Bettys Pond	ACAD	NA	2015-07-30	2717	2015	7	210	SecchiDepth_m	2.25	TRUE	1
ACANTB	Aunt Bettys Pond	ACAD	NA	2018-07-20	3617	2018	7	200	SecchiDepth_m	2.05	TRUE	1
ACANTB	Aunt Bettys Pond	ACAD	NA	2021-07-19	4231	2021	7	199	SecchiDepth_m	1.60	TRUE	1

getSondeInSitu()

This function allows you to query Sonde in situ data by park, site, site type, years, months, parameter, QC type, surface vs. all, etc. Default returns all non-QAQC observations with data. The Sonde_InSitu_Data view is large and can take a few seconds to run for most parks, sites, years, parameters, etc. Surface only measurements, defined as as the medium of all samples within 2m of the surface, are the default.

Get Sonde data for all sites and parameters in MABI from 2021-2023 for non-QAQC samples

mabi <- getSondeInSitu(park = "MABI", years = 2021:2023)
print_head(mabi)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	SiteType	year	month	doy	datetime	QCtype	SampleDepth_m	Parameter	Value	ValueFlag	FlagComments
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	DO_mgL	8.660	NA	NA
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	DOsat_pct	97.200	NA	NA
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	SpCond_uScm	186.400	NA	NA
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	Temp_C	21.000	NA	NA
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	Turbidity_FNU	0.355	NA	NA
MABIPA	The Pogue	MABI	NA	2021-05-24	Lake	2021	5	143	2021-05-24 12:00:00	ENV	0.7445	pH	7.880	NA	NA

get data for all DO parameters in MIMA from 2006-2023 for non-QAQC samples

params <- c("DOsat_pct", "DOsatLoc_pct", "DO_mgL")
period <- 2006:2023
mima_do <- getSondeInSitu(park = "MIMA", years = period, parameter = params)
print_head(mima_do)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	SiteType	year	month	doy	datetime	QCtype	SampleDepth_m	Parameter	Value	ValueFlag	FlagComments
MIMASA	Mill Brook	MIMA	NA	2006-05-15	Stream	2006	5	134	2006-05-15 12:00:00	ENV	0.130	DO_mgL	7.35	NA	NA
MIMASA	Mill Brook	MIMA	NA	2006-05-15	Stream	2006	5	134	2006-05-15 12:00:00	ENV	0.130	DOsat_pct	64.90	NA	NA
MIMASA	Mill Brook	MIMA	NA	2006-06-12	Stream	2006	6	162	2006-06-12 12:00:00	ENV	0.491	DO_mgL	6.38	NA	NA
MIMASA	Mill Brook	MIMA	NA	2006-06-12	Stream	2006	6	162	2006-06-12 12:00:00	ENV	0.491	DOsat_pct	67.40	NA	NA
MIMASA	Mill Brook	MIMA	NA	2006-07-17	Stream	2006	7	197	2006-07-17 12:00:00	ENV	0.231	DO_mgL	5.62	NA	NA
MIMASA	Mill Brook	MIMA	NA	2006-07-17	Stream	2006	7	197	2006-07-17 12:00:00	ENV	0.231	DOsat_pct	64.00	NA	NA

Get Temp data for all sample depths in Jordan Pond in ACAD for non-QAQC samples

ACAD_lake <- getSondeInSitu(site = 'ACJORD', parameter = "Temp_C", 
                            sample_depth = "all")

Get pH for lower NETN parks from May to Oct for QAQC and non-QAQC samples

lnetn_ph <- getSondeInSitu(park = "LNETN", param = "pH", months = 5:10, 
                           QC_type = 'all')

getStreamObs()

This function allows you to query stream observation data by park, site, year and month

Get stream observations for Pogue Stream all years

mabi <- getStreamObs(park = "MABI")
print_head(mabi)

SiteCode	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	WaterClarity	WaterColor	WaterCondition	AlgaeSlimeExtent	AlgaeSlimeAmount	AlgaeSlimeColor	AlgaeFilExtent	AlgaeFilColor	AlgaeClumpsExtent	Algae_Notes
MABIPA	MABI	NA	2013-05-06	2424	2013	5	125			None							NA
MABIPA	MABI	NA	2013-06-04	2423	2013	6	154			None							NA
MABIPA	MABI	NA	2013-07-09	2421	2013	7	189			None							NA
MABISA	MABI	NA	2013-05-06	2136	2013	5	125	Clear	Colorless	None	None			None		None	NA
MABISA	MABI	NA	2013-06-04	2135	2013	6	154	Clear	Colorless	None	None			None		None	NA
MABISA	MABI	NA	2013-07-09	2139	2013	7	189	Clear	Colorless	None	None	Light coating	Greenish	None	Greenish	None	NA

Get observations for all streams in ACAD May 2023

ACAD_streams <- getStreamObs(park = 'ACAD', years = 2023, months = 5)
print_head(ACAD_streams)

SiteCode	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	WaterClarity	WaterColor	WaterCondition	AlgaeSlimeExtent	AlgaeSlimeAmount	AlgaeSlimeColor	AlgaeFilExtent	AlgaeFilColor	AlgaeClumpsExtent	AlgaeClumpsColor	Algae_Notes
ACBRWN	ACAD	NA	2023-05-08	4669	2023	5	127	Clear	Colorless	None	None	N/A	N/A	None	N/A	None	N/A	None
ACCADS	ACAD	NA	2023-05-05	4666	2023	5	124	Clear	Colorless	None	None	N/A	N/A	None	N/A	None	N/A	NA
ACDKLI	ACAD	NA	2023-05-04	4664	2023	5	123	Clear	Colorless	Foamy	None	N/A	N/A	None	N/A	None	N/A	None
ACDUCK	ACAD	NA	2023-05-09	4670	2023	5	128	Clear	Colorless	None	Occasional	Light Coating	Green	None	N/A	None	N/A	None
ACHADB	ACAD	NA	2023-05-05	4667	2023	5	124	Clear	Colorless	None	None	N/A	N/A	None	N/A	None	N/A	None
ACLKWO	ACAD	NA	2023-05-03	4663	2023	5	122	Clear	Colorless	None	None	N/A	N/A	None	N/A	None	N/A	None

getWaterLevel()

This function joins stage and water level data, and allows you to query by park, site, site_type, year, and month. Note that WL data start in 2013 in the view.

Get water level data for Bubble Pond.

bubl <- getWaterLevel(site = "ACBUBL", years = 2013:2023)
print_head(bubl)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	DatumName	DatumType	DatumFunction	Active	TU.TD	StageMethod	DatumLatitude	DatumLongitude	DatumElevation_ft	DatumElevationFeet	GageReadingFeet	WaterLevelFeet	WaterLevel_m
ACBUBL	Bubble Pond	ACAD	NA	2013-05-31	2128	2013	5	150	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.49	332.90	101.4679
ACBUBL	Bubble Pond	ACAD	NA	2013-06-24	2129	2013	6	174	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.77	332.62	101.3826
ACBUBL	Bubble Pond	ACAD	NA	2013-07-18	2130	2013	7	198	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.96	332.43	101.3247
ACBUBL	Bubble Pond	ACAD	NA	2013-08-28	2370	2013	8	239	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.96	332.43	101.3247
ACBUBL	Bubble Pond	ACAD	NA	2013-09-19	2389	2013	9	261	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.78	332.61	101.3795
ACBUBL	Bubble Pond	ACAD	NA	2013-10-24	2410	2013	10	296	SD1	Drilled hole	Stage Measurement	TRUE	TD	Ruler	44.34948	-68.24053	333.39	333.39	-0.78	332.61	101.3795

Get water level data for Weir Pond in August.

weir <- getWaterLevel(site = "WEFAPA", months = 8)
print_head(weir)

SiteCode	SiteName	UnitCode	SubUnitCode	EventDate	EventCode	year	month	doy	DatumName	DatumType	DatumFunction	Active	TU.TD	StageMethod	DatumLatitude	DatumLongitude	DatumElevation_ft	DatumElevationFeet	GageReadingFeet	WaterLevelFeet	WaterLevel_m
WEFAPA	Weir Pond	WEFA	NA	2013-08-28	2226	2013	8	239	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-0.80	568.20	173.1874
WEFAPA	Weir Pond	WEFA	NA	2014-08-12	2668	2014	8	223	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-0.96	568.04	173.1386
WEFAPA	Weir Pond	WEFA	NA	2015-08-20	2827	2015	8	231	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-1.19	567.81	173.0685
WEFAPA	Weir Pond	WEFA	NA	2016-08-25	3161	2016	8	237	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-1.05	567.95	173.1112
WEFAPA	Weir Pond	WEFA	NA	2017-08-24	3403	2017	8	235	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-1.20	567.80	173.0654
WEFAPA	Weir Pond	WEFA	NA	2018-08-14	3664	2018	8	225	SD1	Bolt	Stage Measurement	TRUE	TD	Ruler	41.26022	-73.45142	569	569	-0.51	568.49	173.2758

Plotting Functions

plotWaterBands()

This function produces a plot that summarizes the range of historic data compared with current measurements. The function can handle chemistry, Sonde in situ, Secchi depth, Light Penetration, and Discharge data, although it functions best with sonde and chemistry data. Historic measurements are displayed as min-max values ever previously recorded (outermost band), upper and lower 95% distribution and middle 50% distribution (inner quartiles) of values previously recorded (inner bands). The line represents the median value.

Currently you can only specify one parameter at a time.Values that exceed water quality thresholds (where they exist) are plotted as orange and will show an orange point in the legend. Values within WQ thresholds or for parameters without set thresholds are black. You can include threshold lines (default), or remove them, where they make the y axis range too big, via threshold = FALSE. If multiple sites are specified, they will be faceted.

Plot pH in Jordan Pond for 2023

plotWaterBands(site = "ACJORD", year_curr = 2023, years_historic = 2006:2022, 
  parameter = "pH", legend_position = 'right')

Plot TN in Jordan Pond for 2023, including censored

plotWaterBands(site = "ACJORD", year_curr = 2023, years_historic = 2006:2022,
  parameter = "TN_mgL", legend_position = 'right', include_censored = T)

Same as above, but drop threshold lines

plotWaterBands(site = "ACJORD", year_curr = 2023, years_historic = 2006:2022,
  parameter = "TN_mgL", legend_position = 'right', include_censored = T, threshold = F)

Plot TN in SARA sites in 2023

plotWaterBands(park = "SARA", year_curr = 2023, years_historic = 2006:2022, parameter = "TN_mgL",
  legend_position = 'right')

plotTrend()

This function produces a smoothed trend plot filtered on park, site, year, month, and parameter. It works with chemistry, Sonde in situ, Secchi depth, Light Penetration, and Discharge data. If multiple sites are specified, they will be plotted on the same figure. If multiple parameters are specified, they will be plotted on separate figures.

There are several arguments to customize plots.

Choose whether to include only active sites (default) or all sites that have been monitored via active.
Choose whether to add points, lines, or both (default) via layers argument.
If lines are chosen as a layer, choose whether to plot a loess smoothed line (default) or a line that connects the sample points via smooth.
Choose whether to plot any water quality thresholds that exist via threshold. Upper limits are dashed. Lower limits are dotted.
Choose whether to plot surface (default) or all depth measurements via sample_depth.
Choose whether to include censored values or not via include_censored.
Choose color palette via palette. Default is ‘viridis’, but other options are discrete palettes in RColorBrewer, including ‘Set1’, ‘Dark2’, and ‘Accent’.
Choose position of legend via legend_position. If you don’t want to show the legend, legend_position = 'none'.

Discrete RColorBrewer palettes are below:

RColorBrewer::display.brewer.all(type = 'qual')

Single site; single parameter

Plot non-smoothed surface pH for Eagle Lake for all years.

plotTrend(site = "ACEAGL", parameter = "pH", smooth = FALSE)

Plot smoothed surface pH for Eagle Lake for past 3 years using default span of 0.3 and by default not including the legend.

plotTrend(site = "ACEAGL", parameter = "pH", palette = 'Dark2', years = 2021:2023)

Plot smoothed surface pH for Eagle Lake for all years, removing the legend and using span of 0.75.

plotTrend(site = "ACEAGL", parameter = "pH", span = 0.75)

Plot smoothed Secchi Depth in Jordan Pond for all years, including the legend, different color palette, and using span of 0.75.

plotTrend(site = "ACJORD", parameter = "SDepth_m", 
          span = 0.75, palette = 'Set1')

Plot unsmoothed SO4 in Witch Hole Pond for all years, including censored values.

plotTrend(site = "ACWHOL", parameter = "SO4_ueqL", 
          smooth = F, include_censored = TRUE, legend_position = 'bottom')

Plot unsmoothed line only for pH in Witch Hole Pond for all years.

plotTrend(site = "ACWHOL", parameter = "pH", layers = "lines", 
          smooth = F, legend_position = 'none')

Multiple sites or params

Plot smoothed surface pH for active SARA streams over all years with 0.6 span.

plotTrend(park = "SARA", site = c("SARASA", "SARASC", "SARASD"), 
          site_type = "stream", parameter = "pH", 
          legend_position = "right", span = 0.6)

Plot smoothed surface SO4 for all MIMA streams over all years with 0.6 span.

plotTrend(park = "MIMA", site_type = "stream", 
          parameter = "SO4_ueqL", legend_position = "right", span = 0.6)

Plot non-smoothed surface of multiple Sonde parameters for all MIMA streams over all years with 0.6 span.

params <- c("Temp_C", "SpCond_uScm", "DOsat_pct", "pH")
plotTrend(park = "MIMA", site_type = "stream", 
          parameter = params, legend_position = "right", span = 0.6)

Plot smoothed surface Secchi Depth, Penetration Ratio, pH, and DOC in Jordan Pond for all years, including the legend, different color palette, and using span of 0.75.

plotTrend(site = "ACJORD", parameter = c("SDepth_m", "PenetrationRatio", "pH", "DOC_mgL"), 
          span = 0.75, palette = 'Set1')

Plot smoothed surface water depth in the Pogue for all years, including a different color palette, and using span of 0.75.

plotTrend(site = "MABIPA", parameter = "WaterLevelFeet", span = 0.75, 
          palette = 'Set1')

Plot smoothed TN, TP and SO4 in all MORR sites for all years, including the legend, different color palette, and using span of 0.6. Thresholds, where they exist, plot by default.

plotTrend(park = "MORR", parameter = c("TN_mgL", "TP_ugL", "SO4_ueqL"), 
          span = 0.6, legend_position = 'bottom', palette = 'accent')

plotLakeProfile()

This function produces a heatmap in 1-m bins for ACAD and 0.25m bins for LNETN. You can filter on park, site, year, month, Sonde in situ parameter and either sample relative to the surface or relative to surface elevation. You can only specify one parameter at a time. If multiple sites or years are selected, plots will be faceted on those factors. Keep options limited for best plotting. Note also that you can either select ACAD or LNETN parks, because of the differences in binning.

The option to plot relative to surface elevation (eg depth_type = ‘elev’) corrects sample depth for elevation using water level data and datum elevation for that sampling event. The elevation-corrected option allows you to see how the water column is shifting over time, but currently only works for years >= 2013 in ACAD, as water levels prior to that are not in the current data package. Default setting is raw, where the raw sample depths are plotted instead of elevation.

Note that occasionally profiles skip a bin, which show up as white sections in the plots. Incomplete sampling in 2020 and 2021 also shows up as white sections. If you specify a lake x year x parameter combination that doesn’t exist (e.g., a year a lake isn’t sampled), the function will return an error message instead of an empty plot.

The width of the profiles take into account the number of days between sampling events. For the first and last months (typically May and October), the left/right side of the profiles are padded by 14 days. Otherwise, profile widths are centered on the sample day with the left side representing half the number of days between that visit and the previous visit and the right side representing half the number of days between that visit and the following visit. Black lines are the thermocline, as calculated by rLakeAnalyzer.

There are several arguments to customize plots.

Choose whether to plot the theromocline as points on each profile via plot_thermocline = TRUE (default). The thermocline is calculated by rLakeAnalyzer, and is the depth/elevation at which the largest change in temperature occurs in the sampled water column. If no thermocline is detected, as defined by rLakeAnalyzer::thermo.depth(), nothing is plotted.
Choose whether to include only active sites (default) or all sites that have been monitored via active.
Choose palette. Current enabled themes are ‘viridis’ (yellow - green - blue), and built in continuous color patterns in RColorBrewer. If you prefer other palettes, I can add those too. The only thing I’m trying to avoid is creating the palette manually, since number 1-m bins varies by site and across years. More info on built in ggplot scales can be found here: https://ggplot2-book.org/scales-colour.
Choose position of legend via legend_position. If you don’t want to show the legend, legend_position = 'none'.
Include site name as plot title (title = TRUE). Only enabled when 1 site is selected. Otherwise site names will be in the facets.

Continuous RColorBrewer palettes are below:

RColorBrewer::display.brewer.all(type = 'div')

RColorBrewer::display.brewer.all(type = 'seq')

Simple plots

Plot temperature for Upper Hadlock for years 2013 - 2023 corrected by elevation with thermocline plotted as black lines.

plotLakeProfile(site = "ACUHAD", parameter = "Temp_C", depth_type = 'elev', 
                years = 2013:2023)

Plot temp using raw sample depth (default) for all LNETN lakes sampled in 2023.

plotLakeProfile(park = "LNETN", parameter = "Temp_C", years = 2023, palette = "Spectral")

Plot temperature for Eagle Lake for years 2006 - 2023 with raw sample depth. Note that we can go back to 2006 because we’re using raw sample depth instead of elevation.

plotLakeProfile(site = "ACEAGL", parameter = "Temp_C", depth_type = 'raw', 
                years = 2006:2023)

Same plot as above, but with no plot title or thermocline.

plotLakeProfile(site = "ACEAGL", parameter = "Temp_C", depth_type = 'raw', 
                years = 2006:2023, plot_title = FALSE, plot_thermocline = F)

Plot temperature for all ACAD lakes sampled in 2023 and raw sample depth.

lakes23 <- c("ACBUBL", "ACEAGL", "ACECHO", "ACJORD", "ACLONG", "ACROUN", 
             "ACSEAL", "ACUBRK", "ACUHAD", "ACWHOL")
plotLakeProfile(park = "ACAD", site = lakes23, parameter = "Temp_C", 
                depth_type = 'raw', years = 2023)

Plot DO all ACAD lakes sampled in 2023 and raw sample depth, using reversed RdYlBu palette. Temperature thermocline is also included as black lines.

plotLakeProfile(park = "ACAD", site = lakes23, 
                parameter = "DOsat_pct", depth_type = 'raw', years = 2023, 
                palette = 'RdYlBu', color_rev = TRUE)

Plot specific conductance for Seal Cove Pond from 2013 to 2023 and sample elevation. Uses spectral palette by default and adds site name as title by default. Thermocline is also plotted as default.

plotLakeProfile(site = "ACSEAL", parameter = "SpCond_uScm", 
                depth_type = 'elev', years = 2013:2023)

Plot pH for Jordan Pond from 2013 to 2023 and sample elevation, using reversed spectral palette.

plotLakeProfile(site = "ACJORD", parameter = "pH", 
                depth_type = 'elev', years = 2013:2023, 
                color_rev = TRUE)

Plot pH for all ACAD lakes sampled in 2022 and 2023 and raw sample depth, using reversed RdYlBu palette.

plotLakeProfile(park = "ACAD", site = c("ACBUBL", "ACEAGL", "ACECHO", "ACJORD", "ACLONG", "ACSEAL", 
                                        "ACUHAD", "ACWHOL"), 
                parameter = "pH", depth_type = 'raw', years = 2022:2023, 
                palette = 'RdYlBu', color_rev = TRUE)

Combining plots

Combine plots for temp, DO, pH, and conductance in Bubble Pond for 2023 using the cowplot package.

To minimize typing, I define the parameters I wanted at the beginning. This allows you to adjust the parameters once (i.e., change site), and run through the rest of the code without having to edit it. I also only included the plot title on the first figure, and turned it off for the rest.

The cowplot package must be installed to use this code. Install the package via install.packages('cowplot'). There are other packages to combine plots, including grid and gridExtra, and the function ggarrage() in ggpubr. I tend to start with cowplot, because it’s easy to use and has a great help page. If I really need to customize a plot (like custom spacing for each plot), then I use grid/gridExtra, which allows for more customization, but is a bit harder to work with.

library(cowplot)
sitecode = "ACBUBL"
sitename = getSites(site = sitecode)$SiteName
year = 2023
mon = 5:10
depth = 'elev'
ptitle = F

tplot <- plotLakeProfile(site = sitecode, parameter = "Temp_C", depth_type = depth,
                         years = year, months = mon, plot_title = ptitle)

doplot <- plotLakeProfile(site = sitecode, parameter = "DOsat_pct", depth_type = depth,
                          years = year, months = mon, color_rev = T, plot_title = ptitle)

pHplot <- plotLakeProfile(site = sitecode, parameter = "pH", depth_type = depth,
                          years = year, months = mon, palette = "ryb", color_rev = T,
                          plot_title = ptitle)

cnplot <- plotLakeProfile(site = sitecode, parameter = "SpCond_uScm", depth_type = depth,
                         years = year, months = mon, palette = 'rb',
                         plot_title = ptitle)

# Default settings
plot_grid(tplot, doplot, pHplot, cnplot)

Same plot as above, but customize plot widths, so DO and SpCond have more space for legend, and add title above the grid.

In this case, we’re creating the plot grid with relative widths for each plot. Then we’re combining the title and the plot grid in another grid, setting the title height to be much smaller than the plot grid.

Note that sitename and year are defined in code chunk above.

title <- ggdraw() + draw_label(paste0(sitename, " (", year, ")"), size = 11, fontface = 'bold',
                               x = 0.05, hjust = 0, vjust = 0) 
pgrid <- plot_grid(tplot, doplot, pHplot, cnplot, rel_widths = c(0.9, 0.95, 0.85, 1)) # doesn't seem to be working 
plot_grid(title, pgrid, ncol = 1, rel_heights = c(0.1, 1))

plotPrecipDischarge()

This function produces a plot with dual y-axes, one for precipitation and one for discharge. The x-axis is date. This function only works for one stream monitoring site at a time. Note that ggplot tends to have a lot of warnings that are hard to suppress, particularly for this plot, which has a daily value for precipitation and only monthly values for discharge. You’re always going to get a message about that by using this function. Function currently only plots years where discharge is collected. Discharge is plotted as points and not lines, because discharge can change a lot between sampling events.

Plot Discharge for Mill Brook in MIMA for past 5 years using default colors.

plotPrecipDischarge(site = c("MIMASA"), years = 2019:2023)

Plot Discharge for Aunt Betty Inlet and Kebo Stream for 2022 using different colors. Note that this can be slow because has to download precip. data from NADP. LNETN parks download from a faster web service.

plotPrecipDischarge(site = c("ACABIN", "ACKEBO"), years = 2022, colors = c("cornflowerblue", "orange"))

plotClimTrend()

This function produces a line or smoothed trend plot filtered on park, year, month, and climate parameter. If multiple parks are specified, they can either be plotted on the same figure or separate figures. If multiple parameters are specified, they will be plotted on separate figures. For multiple parks to plot on the same figure, use facet_site = FALSE. For multiple parameters on the same figure, facet_param = FALSE. Note that all months (1:12) are default. If you only want sample months, specify months = 5:10.

NEW: I updated the compiled climate data to cover 1895 up through May 2024, so you can plot historic trends beyond just the 2006-2024 monitoring period.

Plot avg temp for LNETN from 1895:2024, without points and with smoothed line and span 0.1.

plotClimTrend(park = "LNETN", years = 1895:2024, layers = 'lines', 
              parameter = "tmean", span = 0.1, palette = "Dark2", legend_position = 'right')

Plot monthly precip for MORR, ROVA and WEFA from 1895:2024, without points and with smoothed line and span 0.1.

plotClimTrend(park = c("MORR", "ROVA", "WEFA"), years = 1895:2024, layers = 'lines', 
              parameter = "tmean", span = 0.1, palette = "Dark2", legend_position = 'right')

Plot temp stats from 1895-2024 for ACAD on the same figure.

plotClimTrend(park = "ACAD", years = 1895:2024, parameter = c("tmin", "tmean", "tmax"), span = 0.1, layers = 'lines', facet_param = F, legend_position = 'right', palette = "Set2")

Plot monthly climate data for MABI from 2006:2023 (all months is default), without a smoothed line.

plotClimTrend(park = "MABI", years = 2006:2023, parameter = "all", smooth = F, facet_param = T)

Plot same monthly NOAA data for MABI from 2006:2023, with smoothed line and span 0.7, and using the Dark2 color palette.

plotClimTrend(park = "MABI", years = 2006:2023, parameter = "all", span = 0.7, palette = "Dark2")

Plot monthly NOAA min, max and mean temperature for MABI and SARA from 2006:2024 (May), with smoothed line, span 0.7, and parks on separate plots, but parameters on the same plots.

plotClimTrend(park = c("MABI", "SARA"), years = 2006:2023, 
              parameter = c("tmean", "tmax", "tmin"), 
              span = 0.7, facet_park = TRUE, facet_param = FALSE,
              legend_position = "bottom")

plotClimAnom()

This function produces a trend plot filtered on park, year, month, and climate parameter and shows anomalies from chosen baseline (1901 - 2000 or 1991 - 2020) as red for above and blue for below average using gridded NOAA NClim data. If a year x month combination is specified that doesn’t occur yet in NETN_clim_annual dataset, it will be downloaded if available. New months are typically available within a few weeks of the month end. If multiple parks are specified, they can either be plotted on the same figure or separate figures. If multiple parameters are specified, they will be plotted on separate figures.

Plot ACAD average temperature anomalies from 1895 to 2024 compared with 20th century normals.

plotClimAnom(park = "ACAD", years = 1895:2024, parameter = "tmean", legend_position = 'right')

Plot MABI and SAGA average temperature anomalies from 2006 to 2024 compared with 20th century normals.

plotClimAnom(park = c("MABI", "SAGA"), years = 2006:2024, parameter = "tmean", legend_position = 'right')

Plot MABI and SAGA precip. anomalies from 2023 compared with 20th century normals.

plotClimAnom(park = c("MABI", "SAGA"), years = 2023, parameter = "ppt", legend_position = 'right')

Plot ACAD avg temp and precip anomalies from last 20 years compared with 30-year normals.

plotClimAnom(park = c("ACAD"), years = 2004:2024, parameter = c("tmean", "ppt"), 
             legend_position = 'right', averages = "norm1990")

plotClimComps()

This function compares either the 20th century normals (1901 - 2000), or the most recent 30-year normal (1991 - 2020) average monthly climate variables with user-specified years to provide an idea of how extreme or normal a given month in a year is.

This function has been overhauled to work with NOAA gridded climate data for both the normals and the current data. The same algorithms are used to generate historic and current data, which minimizes the bias introduced by weather station vs. Daymet data (Prism had similar issues). Using the NOAA dataset is slower to process, so I went with park-level centroids to extract climate data instead of extracting at the site level. I also compiled all of the NOAA data from 1895 through May 2024. I have a webalert that tells me when a new month is available, so I can keep the compiled data up to date (see data("NETN_clim_annual")). However, if I get behind, this function will download the missing months of data that are available online.

Single parameter

Plot mean monthly temp for MABI for 2019:2023 and all months with red-blue color palette and the 4 decadal averages.

plotClimComps(park = "MABI", years = 2019:2023, parameter = "tmean", 
              palette = c("red", "blue"), averages = "norm20cent")

Same as above, but with points and 30-year normal starting at 1991.

plotClimComps(park = "MABI", years = 2019:2023, parameter = "tmean", palette = c("red", "blue"), 
              layers = c('points', 'lines'), averages = "norm1990")

Plot max monthly temp for ROVA for 2023 and 2024, with ‘viridis’ palette and 20th century norm (default).

plotClimComps(park = "ROVA", years = 2023:2024, parameter = "tmax", 
              palette = "viridis", layers = 'lines')

Plot total monthly precip for Jordan Pond for past 5 years using orange-blue color palette, compared with the 30-year norm starting at 1990.

plotClimComps(park = "ACAD", years = 2018:2023, parameter = 'ppt', 
              palette = c("orange","blue"), averages = "norm1990")

Plot total monthly precip for ACAD for current year using red for 2023 and blue for 2024

plotClimComps(park = "ACAD", years = 2023:2024, parameter = 'ppt', palette = c("red", "blue"))

Plot total monthly max temp for ACAD for current year using red for current year color and 30-year normal.

plotClimComps(park = "ACAD", years = 2024, parameter = 'tmax', palette = c("red"), averages = "norm1990")

Multiple plots

The code below would take a lot to program into waterNETN functions, so I’m instead showing you how you can create panels of multiple plots using the cowplot package.

Generate grid of all parameters in ACAD from past 5 years

Note that the code that names leg is extracting a legend from a copy of the figures to be plotted separately from the other figures. This allows you to control the sizing and alignment better. The NULLs also allow you to create white space between plots.

tmax <- plotClimComps(park = "ACAD", years = 2019:2023, parameter = "tmax",
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"))

tmin <- plotClimComps(park = "ACAD", years = 2019:2023, parameter = "tmin", 
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"))
tmean <- plotClimComps(park = "ACAD", years = 2019:2023, parameter = "tmean",  
                       layers = 'lines', legend_position = 'none', palette = c("red", "blue"))
prec <- plotClimComps(park = "ACAD", years = 2019:2023, parameter = 'ppt', 
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"))

prec_leg <- plotClimComps(park = "ACAD", years = 2019:2023, parameter = 'tmax', 
                          layers = 'lines', legend_position = 'right', palette = c("red", "blue")) +
 theme(legend.box.margin = margin(0, 0, 0, 5), 
       legend.text = element_text(size = 9), 
       legend.title = element_text(size = 9))

leg <- cowplot::get_plot_component(
  prec_leg,
  'guide-box-right',
  return_all = T
  )


grid1 <- cowplot::plot_grid(tmean, tmax, NULL, tmin, prec, NULL, 
                            nrow = 2, ncol = 3, rel_widths = c(1, 1, 0.05, 1, 1, 0.05))
grid2 <- cowplot::plot_grid(grid1, NULL, leg, ncol = 3, nrow = 1, rel_widths = c(1, 0.05, 0.3))
grid2

Generate grid of all parameters in ACAD for 2024, Jan through May.

Because waterNETN plotting functions return ggplot objects, we can tweak the plots so that temp plots all have the same y axis (so you can see that max is higher than min).

data("NETN_clim_annual") # compiled annual data
data("NETN_clim_norms") # compiled normals
head(NETN_clim_annual)

##   UnitCode
## 1     ACAD
## 2     BOHA
## 3     ROVA
## 4     MABI
## 5     MIMA
## 6     MORR
##                                              UnitName
## 1                                Acadia National Park
## 2      Boston Harbor Islands National Recreation Area
## 3        Roosevelt-Vanderbilt National Historic Sites
## 4 Marsh-Billings-Rockefeller National Historical Park
## 5                 Minute Man National Historical Park
## 6                 Morristown National Historical Park
##        long
## 1 -68.26019
## 2 -71.02159
## 3 -73.92672
## 4 -72.53824
## 5 -71.29604
## 6 -74.53749
##        lat
## 1 44.37656
## 2 42.27097
## 3 41.76804
## 4 43.63382
## 5 42.45340
## 6 40.76780
##        prcp
## 1 150.50000
## 2  96.19043
## 3  92.76953
## 4  55.83984
## 5  92.24023
## 6 118.08984
##         tmax
## 1 -0.5898438
## 2  0.7099609
## 3 -1.7099609
## 4 -2.6503906
## 5  0.1298828
## 6  0.2197266
##         tmin
## 1 -12.139648
## 2  -8.129883
## 3 -10.780273
## 4 -16.750000
## 5 -10.809570
## 6  -9.500000
##        tavg
## 1 -6.370117
## 2 -3.709961
## 3 -6.250000
## 4 -9.700195
## 5 -5.339844
## 6 -4.639648
##   year
## 1 1895
## 2 1895
## 3 1895
## 4 1895
## 5 1895
## 6 1895
##   month
## 1     1
## 2     1
## 3     1
## 4     1
## 5     1
## 6     1

ymin <- round(min(NETN_clim_annual$tmin, NETN_clim_norms$tmin_norm_1901_2000,     
                  NETN_clim_norms$tmin_norm_1991_2020), 0) - 1
ymax <- round(max(NETN_clim_annual$tmax, NETN_clim_norms$tmax_norm_1901_2000, 
                  NETN_clim_norms$tmax_norm_1991_2020), 0) -1

tmax <- plotClimComps(park = "ACAD", years = 2024, parameter = "tmax", data_type = "wstn",
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"),
                      plot_title = F) +
        ylim(ymin, ymax)
tmin <- plotClimComps(park = "ACAD", years = 2024, parameter = "tmin", data_type = "wstn", 
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"),
                      plot_title = F) +
        ylim(ymin, ymax)
tmean <- plotClimComps(park = "ACAD", years = 2024, parameter = "tmean", data_type = "wstn",
                       layers = 'lines', legend_position = 'none', palette = c("red", "blue"),
                      plot_title = F)+
        ylim(ymin, ymax) 
prec <- plotClimComps(park = "ACAD", years = 2024, parameter = 'ppt', data_type = "wstn",
                      layers = 'lines', legend_position = 'none', palette = c("red", "blue"),
                      plot_title = F)

prec_leg <- plotClimComps(park = "ACAD", years = 2024, parameter = 'tmax', data_type = "wstn",
                          layers = 'lines', legend_position = 'right', 
                          palette = c("red", "blue")) +
 theme(legend.box.margin = margin(0, 0, 0, 5), 
       legend.text = element_text(size = 9), 
       legend.title = element_text(size = 9))

leg <- cowplot::get_plot_component(
  prec_leg,
  'guide-box-right',
  return_all = T
  )


grid1 <- cowplot::plot_grid(tmean, tmax, NULL, tmin, prec, NULL, 
                            nrow = 2, ncol = 3, rel_widths = c(1, 1, 0.05, 1, 1, 0.05))
grid2 <- cowplot::plot_grid(grid1, NULL, leg, ncol = 3, nrow = 1, rel_widths = c(1, 0.05, 0.3))
grid2

plotClimCumPrecip()

This function plots cumulative monthly precipitation for a given year compared to either the 19th century or 30-year normal. If multiple parks or years are specified, resulting plots will be faceted on those variables.

Plot ACAD cumulative precipitation for 2020 through 2023

plotClimCumPrecip(park = "ACAD", years = 2020:2023, legend_position = 'bottom')

Plot ACAD cumulative precipitation for 2013 through 2023 for April - October only

plotClimCumPrecip(park = "ACAD", years = 2014:2024, legend_position = 'bottom', months = 4:10)

Plot all but SAIR cumulative precipitation for 2023 with 4 columns

parks <- c("ACAD", "MABI", "MIMA", "MORR", "ROVA", "SAGA", "SARA", "WEFA")
plotClimCumPrecip(park = parks, years = 2023, legend_position = 'bottom', numcol = 4)

Plot cumulative precipitation for 2019 - 2023

plotClimCumPrecip(park = "ACAD", years = 2019:2023, legend_position = 'bottom', numcol = 5)

plotClimRel()

This function plots a given year compared against either the 19th century normal (1901 - 2000) or the 1991 - 2020 normal. Both normals and annual comparisons are derived from NOAA gridded climate data. You can either plot all temp variables, min/max temp, mean temp or precipitation. If multiple years are specified, resulting plots will facet on year. You can make your own palette by specifying colors, or by default use viridis. Currently can only plot one park at a time.

Plot all temperature variables on 1 graph for MABI in 2024.

plotClimRel(park = "MABI", years = 2024, parameter = "temp", palette = c("#EEE55A", "#D56062", "#067BC2"))

Plot precip for SARA 2023 compared to 1990 - 2019 normals.

plotClimRel(park = "SARA", years = 2023, parameter = "ppt", palette = "grey", averages = "norm1990")

Plot precip for SARA 2023 and 2024 compared to 1990 - 2019 normals.

plotClimRel(park = "SARA", years = 2023:2024, parameter = "ppt", 
            palette = "grey", averages = "norm1990")

Plot precip for LNETN parks in 2023 without a legend and with 4 columns instead of 3.

plotClimRel(park = "LNETN", years = 2023, parameter = "ppt", palette = "grey", legend_position = 'none', numcol = 4)

Same as above, but with full park names as facet titles.

plotClimRel(park = "LNETN", years = 2023, parameter = "ppt", palette = "grey", legend_position = 'none', numcol = 4, title_type = "UnitName")

plotClimDrought()

This function plots weekly drought index at the county level for each specified park or weather station nearest to a specified park. Resulting plot shows the percent of county area in 5 levels of drought, with D0 = Abnormally Dry, D1 = Moderate Drought, D2 = Severe Drought, D3 = Extreme Drought, and D4 = Exceptional Drought. Drought designations come from the U.S. Drought Monitor. If multiple parks are specified, results will be faceted with a separate plot for each park. Where multiple counties occur in a park, data can be faceted by county if dom_county = FALSE. To only plot predominant county, specify dom_county = TRUE (default).

I tweaked the x-axis formatting a bit based on whether 1 year or many years are plotted. Let me know if you’d like further tweaks.

Plot drought info for ACAD in 2023

plotClimDrought(park = c("ACAD"), years = 2023)

Plot drought info for MABI and SAGA for 2020 to 2024

plotClimDrought(park = c("MABI", "SAGA"), years = 2020:2024)

Plot drought info for both MORR counties for 2023

plotClimDrought(park = "MORR", years = 2023, dom_county = FALSE)

Plot drought info for both ACAD counties in 2020 for May through October

plotClimDrought(park = "ACAD", years = 2020, dom_county = FALSE, months = 5:10)

Plot drought info for LNETN in 2023 with legend on bottom

plotClimDrought(park = "LNETN", years = 2023, legend_position = 'bottom')

Summary Functions

sumClimAvgs()

Now that we’re using NOAA gridded climate data, we no longer need this function. It’s been removed from the package.

sumClimMonthly()

Now that we’re using NOAA gridded climate data, we no longer need this function. It’s been removed from the package.